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We obtain the exact time evolution for the one-dimensional integrable fermionic 1/r Hubbard 
model after a sudden change of its interaction parameter, starting from either a metallic or a Mott- 
insulating eigenstate. In all cases the system relaxes to a new steady state, showing that the presence 
of the Mott gap does not inhibit relaxation. The properties of the final state are described by a 
generalized Gibbs ensemble. We discuss under which conditions such ensembles provide the correct 
statistical description of isolated integrable systems in general. We find that generalized Gibbs 
ensembles do predict the properties of the steady state correctly, provided that the observables or 
initial states are sufficiently uncorrelated in terms of the constants of motion. 
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I. INTRODUCTION 

Recent experiments with ultracold atomic gases P, d, 
H, [1] have made it possible to study the time evolution 
of a tunable quantum many-body system that is kept in 
excellent isolation from the environment. For example, 
such a quantum system can be forced out of equilibrium 
by suddenly changing a parameter in the Hamiltonian. 
Then the system may or may not relax to a new steady 
state, which is not necessarily the thermal state predicted 
by statistical mechanics. After such a "quantum quench" 
the system evolves according to Schrodinger's equation 



= exp(-ii7i/ft)|*(0)) 



(1) 



where |^(0)) is the prepared initial state and H is the 
new Hamiltonian for times t > 0. This situation has 
recently been studied by a variety of numerical and ana- 
lytical techniques 
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Due to the unitary time evolution the wave function 
l^'(t)) of an isolated system remains pure for all times 
and docs not converge for t — > oo. Only the state of a 
finite subsystem, for which the rest of the system effec- 
tively acts as a reservoir, can become stationary (l^[l6j. 
Nevertheless, also for the entire system we expect relax- 
ation of the expectation value {"^^ {t)\0\^ {t)) of an ob- 
servable O to a stationary value for large times. How- 
ever, this global relaxation can happen only for (i) suffi- 
ciently large systems, (ii) sufficiently simple observables, 
and (iii) sufficiently complicated Hamiltonians, for the 
following reasons. First of all, (i) many degrees of free- 
dom are needed, so that the thermodynamic limit may be 
takeii, otherwise one expects finite recurrence times [T9| 
(see [201 for a recent example). Furthermore, (ii) the 
expectation value of a complicated observable need not 
relax; for example, the expectation value of O = |rii)(ni| 
+ 1^2) {n2\, involving the projectors onto two cigcnstatcs 



of H with different energies, oscillates for all times. Usu- 
ally such projectors are highly nonlocal and their expec- 
tation values correspond to correlation functions of very 
high order. On the other hand, local and few-particle 
observables arc usually simple enough to relax to new 
stationary values. Finally, (iii) the Hamiltonian H that 
governs the dynamics must also be sufficiently compli- 
cated. For example, the magnetization of an Ising chain 
in a transverse magnetic field relaxes for long-range inter- 
actions, but keeps oscillating when only next neighbors 
are coupled [2l[. Similar "coUapse-and-revival" oscilla- 
tions of a many-body system were recently observed in 
experiments with ultracold atoms by Greiner et al. 01- 
In their experiments, a Bosc condensate was prepared 
in the potential of an optical lattice which was suddenly 
steepened 0, then the bosons are essentially only subject 
to the Hubbard interaction H = t^X^i "-iC^-j ~ 1) but no 
hopping between lattice sites occurs. Since in this case 
H /U has only integer eigenvalues, the wave function will 
oscillate for all times with period 2tiU /h. 

For small hopping between lattice sites it follows from 
perturbation theory that expectation values keep oscil- 
lating for short times. But it is not clear what will hap- 
pen for long observation times. Is relaxation possible 
for a bosonic or fermionic Hubbard model if the spec- 
trum has a Mott gap, or if the initial state is a Mott 
insulator? Or is it prevented by the Mott gap in the 
energy spectrum? The answer is no, not necessarily: in 
Sec. im we provide an example, the 1/r fermionic Hubbard 
model [2^, which shows that relaxation in the presence 
of a Mott gap is indeed possible. Note that the formation 
of a fermionic Mott insulator was recently observed with 
ultracold atoms p3l |. 

Another central question is whether the steady state of 
a quenched isolated system can be described by an effec- 
tive density matrix p, such that Tr[Op] yields the correct 
expectation value for any observable O which relaxes. 
Statistical mechanics can be used to make an approxi- 
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mate but usually accurate prediction pmic for this steady- 
state density matrix. For example, the microcanonical 
prediction is that pmic ~ const for states with energy 
close to (^'(0)|i7|*(0)), and zero otherwise. If Tr[Op,„ic] 
indeed agrees with the long-time limit of {'i> {t)\0\'i> (t)) , 
we say that the system thermalizes. Clearly, thermaliza- 
tion can be expected only for sufficiently coarse-grained 
observables; it is always possible to construct a compli- 
cated correlation function that depends on the details of 
the initial conditions and is not described by Pmic- As for 
classical gases, thermalization is generally expected for 
isolated interacting quantum systems. Indeed, for one- 
dimensional atomic Bose gases, the dynamics leading to 
the thermal state were recently observed by Hofferberth 
et al. 

Thermalization of an isolated system is impossible in 
certain cases, usually because the system is integrable in 
the sense that there are infinitely many constants of mo- 
tion. Simple theoretical examples are one-dimensional in- 
tegrable models, such as the XY spin chain, for which the 
magnetization does not thermalize after a quench of the 
longitudinal magnetic field [H, [H, [H, [13, or the Ising 
chain in a transverse field, whose correlation functions 
do not thermalize after a field quench Q. Experimen- 
tally, the lack of thermalization was recently observed 
by Kinoshita et al. [H for bosonic atoms confined to one 
dimension, whose nonthermal stationary momentum dis- 
tributions were attributed to integrability. 

The fundamental assumption of statistical mechanics 
is that the equilibrium state is characterized only by a 
few thermodynamic variables such as internal energy and 
particle number. However if the system is integrable, its 
infinitely many constants of motion lead to a rather de- 
tailed memory of the initial state, because much fewer 
states are accessible to the dynamics. Nevertheless, even 
in the absence of thermalization, a statistical prediction 
for the steady state can be made with a generalized Gibbs 
ensemble ( GGE), as discussed in the recent work of Rigol 
et al. 0. For example, if the time evolution of an inte- 
grable system is determined by the effective Hamiltonian 



(2) 



where the operators Xq commute, [Xq;,X^] = 0, then the 
standard choice for the statistical operator of the GGE 
is constructed from these constants of motion according 
to i 



Pg 



= Tr[e 



(3) 



This choice maximizes the entropy (S = — Tr[pln/3]) for 
fixed expectation values (Xq,)q, which are set to their 
initial-state expectation values (2q,)o by an appropriate 
choice for the Lagrange multipliers Aq [11]. The GGE 
prediction for the steady-state expectation value of an 
observable O is then (0)q = Tr[OpQ]. GGEs success- 
fully predict some properties of the nonthermal steady 



states occurring after quenches in integrable or highly 
constrained systems. For example, they yield the correct 
nonthermal momentum distribution of one-dimensional 
hard-core bosons [^ (experimentally realized in Ref. [^, 
for the one-dimensional Tomonaga-Luttinger model |8|, 
and for the Falicov-Kimball model in infinite dimen- 
sions [l^. The GGE also yields the correct double oc- 
cupation for the 1/r fermionic Hubbard model, as shown 
in Sec. [Til below. 

However, in some cases, the stationary values of some 
observables differ from the statistical predictions of the 
GGE. For one-dimensional hard-core bosons the unit- 
cell averaged one-particle correlation function is not de- 
scribed by the GGE [^. Furthermore Gangardt and 
Pustilnik [l5[ pointed out that the GGE ^ may not 
capture correlations between the conserved quantities ■ 
As a consequence the merit of generalized Gibbs ensem- 
bles is currently somewhat controversial. In this situation 
rather general criteria for the validity of GGEs should be 
useful, which we derive in Sec. IIIII Our approach is com- 
plementary to the work of Barthel and SchoUwock [l6l |. 
who recently showed that for finite subsystems the re- 
duced density matrix converges to the GGE under certain 
mathematical conditions on the initial state and Hamil- 
tonian. 

It should be noted that an important ambiguity lingers 
in the construction of the GGE for the Hamiltonian 
While all the are conserved, so are all their combina- 
tions, i.e., all products of the form Xq,X^, laTpX^, . . . , 
leading to the question of whether all such products 
should also be included in the exponent of the density 
matrix ([3|). In Sec. IIIII wc show that for observables and 
initial states which involve little or no correlation, it suf- 
fices to fix the constraints Xq,, as in Eq. ([3]), but not their 
products. We also provide an example for which different 
choices of the GGE lead to different predictions. 

Finally, we note that thermalization might also be 
prevented in nonintegrable systems due to many-body 
effects, e.g., the presence of a Mott gap. Nonther- 
mal steady states in nonintegrable systems were ob- 
served and argued for in recent numerical studies for fi- 
nite one-dimensional soft-core bosons [l3| and spinless 
fermions [llj. By contrast, thermalization was observed 
for hard-core bosons on a two-dimensional lattice [l^ . 
Fast relaxation to a nonthermal quasisteady state, so- 
called prethermalization [2^ . was recently observed for 
the fermionic Hubbard model in high dimensions [l7j . 
Further studies of relaxation in nonintegrable many-body 
systems are therefore desirable, but will not be the sub- 
ject of this paper. 

Our goals in the present paper are thus two-fold, (i) 
On the one hand, we provide an explicit example of relax- 
ation in a fermionic Mott insulator: In Sec.[TT|we obtain 
the exact time evolution of the one-dimensional fermionic 
Hubbard model with 1/r hopping and repulsive interac- 
tion U [13 , starting from the metallic ground state at U 
= or the insulating ground state at [/ = oo. Wc find 
that the expectation value of the double occupation d{t) 



3 



relaxes with algebraically damped oscillations to a new 
stationary value for all U, i.e., relaxation is not inhibited 
by the presence of a Mott gap. The long-time limit doo = 
limi_,oo d{t) differs from the thermal value, as expected 
for an integrable system, but is described by an appropri- 
ate GGE. (ii) On the other hand we discuss under which 
circumstances GGEs describe the steady state of inte- 
grable systems after relaxation. We show in Sec. IIIII that 
their validity depends on the observable, on correlations 
in the initial state, and possibly the system size. 



II. INTERACTION QUENCH IN A FERMIONIC 
HUBBARD MODEL 

1/r Hubbard chain 

We consider sudden changes in the interaction parame- 
ter of the one-dimensional 1/r fermionic Hubbard model. 



Ho 



(4a) 
(4b) 



ij = l..L 



|fc|<7r 



with repulsive on-site interaction U, bandwidth W, and 
dispersion ek=Wk/{2'K), which corresponds to hopping 
amplitudes t„y = {-iW/2L){-iy"~^ / sm[7:{m - j)/L] 
that decay proportionally to inverse distance. This lat- 
tice model was introduced by Gebhard and Rucken- 
stein [13] as a parent system of the 1 /r^ Haldanc-Shastry 
Heisenberg chain [s^, [3l|, to which it reduces in the 
limit of large U for a half-filled band with density n 
= 1. We consider only n < 1; larger densities can be 
treated by means of a particle-hole transformation, Cf.^ 

~* (^-ktr Si- For U > —W and any number of lattice 
sites L, the model ^ is represented by an effective nonin- 
teracting bosonic Hamiltonian, from which ground-state 
and thermodynamic properties can be obtained analyti- 
cally [13, m, . For [/ = the ground state of ^ is of 
course the Fermi sea, 



iv-o) = n 



h<CkF 



with particle density n = 1 + k-pjj:. For U 
other hand, the ground state is |2^ . [33 | 



(5) 



oo, on the 



(6) 



i.e., the Fermi sea with all doubly occupied sites projected 
out. At half-filling a Mott-Hubbard metal- insulator tran- 
sition occurs at interaction strength Uc = W , with the 
Mott gap given by A = [/ - C/c for [/ > Uc [H]- This 
metal-insulator transition is also captured by correlated 
variational wave functions [3J, Sal ■ 



Interaction quenches 

We now consider the following nonequilibrium situa- 
tion. For times t <Q the system is prepared in the ground 
state for interaction parameter U = or oo, i.e.. 



I*(0)) - 
l*(0)> - 



IV'oo) 
v/(V'oo|V'oo) 



metallic state, or (7a) 
Mott insulator. (7b) 



Then at time i = the interaction is suddenly switched 
to a new value < J7 < oo, so that the time evolu- 
tion for < > is governed by the Hamiltonian (|3]), i.e., 
the system evolves according to Eq. ([ij. We refer to 
these two types of quenches as ^ C/ (starting from the 
metallic state ((7a|)) and oo — > J7 (starting from the Mott- 
insulating state (j7bp ). respectively. More general initial 
states corresponding to intermediate values of U can also 
be used; they lead to similar results which are omitted 
here. 



Bosonic representation 

In the bosonic representation of Ref. [sj the initial 
states ([7]) factorize. They can be written in terms of 
hard-core bosons (•, o, "f , |) in the form [13, [s^ 



\%) = I [Ti]---[Ti] 

K 



(8a) 
(8b) 

(8c) 



where A = 27r/i, /Cp = (2n — l)7r, and the prime in- 
dicates that only every other bosonic pseudomomentum 
K, G {—T^, tt) appears. The Hamiltonian ^ only acts 
separately on each space spanned by the bracketed con- 



figurations [t I] = (p) and 



pseudomomenta, K. and K. + A [22, l34 



for neighboring 
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W ( sgn(/C) 







-sgn(/C) 



D 



1 - a/C \n-o?K 



(9a) 

(9b) 
(9c) 



with ajc = sgn(/C)(2/C + A)/{2tt). In this representa- 
tion it is then straightforward, although tedious, to ob- 
tain the propagator exp{—iHt), its action on |^'(0)), and 
the expectation value of the double occupation d{t) = 
(^(t)|nif nij^(t)). We take the thermodynamic limit, 
L — > oo with fixed density n; the sums over /C are then 
replaced by integrals which can be evaluated analytically. 
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Results for the double occupation 

Setting the bandwidth to W = 1 (and also h ~ 1), our 
results for the interaction quenches ^ U and oo s- U 
can be written as 



dit) 
d{t) 



c+ + fit) 
c--fit) 



for the two types of quenches, with the abbreviations 

,2 ^2 



(10a) 
(10b) 



n 



32t/2 



2nU + In - 



n Lot sin(a;t) + 3 coslujt) 
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(11) 
(12) 



These expressions involve several energy scales, apart 
from the interaction U and the bandwidth W (= 1), 
namely the Mott gap A = U — 1, the total bandwidth of 
the spectrum il = 2 + A = U+l, and uj = ^/fl'^ - iUn, 
a characteristic density-dependent energy scale appear- 
ing in the holon and spinon excitation energies |22l . [ssj . 
As functions of U the constants c± have the remarkable 
symmetry that both are invariant under the replacement 
U 1/U (for all n). The function g{r],t) in Eq. ^ is 
given by 



1 



32C/2 



t3 



7/sin(77t) 



<4 



cos{r]t) , (13) 



where Ci(a:;) = — JJj" cos{y)/y is the integral cosine. For 
the quench to U = Uc (= 1) this reduces to c± = 1/8 
and f{t) = fi{t), where 



2t' 
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3 . ^ 

^sm(2t)- 



16*4 



■cos(2i). (14) 



We note that in all cases d{t) relaxes with damped 
oscillations from its initial value d(0) (= n'^/A for the 
metallic state ((7a)) or for the Mott insulator (|7bp ) to 
a new stationary value. This long-time limit, d^ = 
Ihnt— too d(t) , always exists, even when quenching to U 
= 0, [/c, or oo. For these final values of U we find 



lim d{t) 

t — *oo 



n 
T 



n2(3-2n) 
6 



U 
U 

u 











(15) 



1 or cx) 



For the quench to U = we note in particular that the 
stationary value of the double occupation differs from the 
thermal value, -n? /A. This is discussed in detail at the end 

of Sec. nni 



For general U and n the function (|12p behaves asymp- 
totically as 



fit) 



ri(l — n) sin(cj*) cos(r2<) 



LOt 



(1 - 3n)cj2 + (1 



4f/i2 
n)ri2 cos(a;t) 



2^2 



4f/t2 



(16) 



for large times, which at half-filling reduces to perfect 
beating, 



1 

8^ 



(l-i7)2 (l-[/2)2 



16C/2 



In 



1 - U 



l + U 



fit) - - 



cos(Ot) + cos(Ai) cos([/t) cos(t) 



4[/<2 



2[/i2 



(17) 
(18) 



We conclude that the relaxation of dit) involves the fre- 
quencies to and Q and that it falls off rather slowly as 
1/t for densities n < 1, or as 1/^2 for n = 1. This 
type of algebraic decay is typical for one-dimensional 
systems [33 • Mathematically it can be traced to the 
Riemann-Lebesgue lemma, i.e., the 0(l/<) decay of one- 
dimensional integrals over oscillating functions 



Nonthermal steady states 
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The results for dit) for several quenches are shown 
in Fig. [1] together with the long-time limit (dotted 
blue lines) and the thermodynamic prediction (solid red 
lines). The latter is determined from the exact grand- 
canonical potential f{T,fi,U) = -iT/L) Tr[exp(-(i/ - 
^N)/T)] m,!!!,!!!], using the temperature T and chem- 
ical potential fi that correspond to the same internal en- 
ergy and density as the final state. The density is given 
by {N)gca,n/L = df/dfi, which yields the chemical po- 
tential /i(T, n, U) by inversion, and the internal energy 
per site is eiT,^i,U) = (if)gca„/i = -T^dif/T)/dT - 
df /dfi. The temperature T is obtained as the solution of 
e(r, ^(T, n, U), U) ~ {H)t>o, where the energy in the fi- 
nal state is {H)t>o = — (2 — n)7i/4 + C/2/4 for the metallic 
state (|7ap and {H)t>o = —(1 — n)n/4 for the insulating 
state (|7bp . Finally, the thermal value of the double occu- 
pation per site is obtained as (i(T, /i, U) — ('^j|'^jx)gcan = 
df /dU, evaluated at the determined temperature T and 
chemical potential /x(r, n, U). 

The long-time limit of dit) clearly differs from the ther- 
modynamic prediction (see Fig.[T]). This is the expected 
behavior for an integrable system with an infinite number 
of constants of motion. They act as constraints on the 
accessible states and must be taken into account into the 
statistical description of the steady state. As discussed 
in the Introduction, this can be achieved by employing 
a GGE For the 1/r Hubbard chain (g]) the effective 
bosonic Hamiltonian is indeed of the form where a 
labels the pairs of bosonic pseudomomenta (/C,/C -I- A) 
and the constants of motion Xa arc the projectors onto 
one of the two eigenstates \i'icx+^) of ^tc = T/c^+Aic 
+ UDk.,k+ak [Eq. [5]. For the initial states © of the 
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FIG. 1: (Color online) Relaxation of the double occupation d{t) after an interaction quench to U — 0.5 (top half) or 1.5 
(bottom half), for density n = 1 (left column) and n = 0.7 (right column), from the metallic state (|7a|l in (a,c,e,g) and from 
the insulating state H7b|) in (b,d,f,h). The dashed black line shows the leading-order envelope of the asymptotic behavior. The 
leading-order asymptote is marked by the dashed green line and the long-time limit by the dotted blue line. The solid red line 
marks the grand-canonical prediction for a temperature and chemical potential that corresponds to the internal energy and 
density of the final state. The insets shows the large-t behavior in more detail. 
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form l^'(O)) = Yi'fc IV'K./c+a) the statistical operator of 
the GGE (El) becomes 



Pg = IT' X] I(^'c,/c+a|^/'k;,/c+a)P Wk,k.+a){i^k,ic+a\ , 

(19) 



K iy=l 



from which we obtain {niinii)Q as c+ and C-/U for the 
two types of quenches. The double occupation that is 
predicted by the GGE thus agrees precisely with the long- 
time limit [Eq. i.e., {niiriii)^ = liiat^ao d(t). In 
the next section we discuss general rules when GGEs are 
valid, which will explain in particular why the GGE gives 
the correct stationary value of the double occupation in 
the 1/r Hubbard model. 



III. VALIDITY OF GENERALIZED GIBBS 
ENSEMBLES 

Long-time average and diagonal ensemble 

In order to address the validity of GGEs we first need 
to determine the long-time limit of a quantum system 
after an arbitrary quench. Suppose that an isolated sys- 
tem is prepared at time < = in an initial state which 
is described by the density matrix po, while the time 
evolution for t > is governed by an arbitrary time- 
independent Hamiltonian H. For an initial pure state 
(e.g., as in Eq. ([7])) the density matrix po has the form 
Po = |*(0))(«'(0)|, whereas po = E„Pn|*n)(*n| for a 
statistical mixture of orthogonal states |4'„) with proba- 
bilities Pn (e.g., for an initial grand-canonical ensemble, 
as in Ref . [l3l ) . 

For t > the time evolution of the density matrix is 
given by 



Pit) = e^^Voe- 
and the expectation value of an observable O is 



(20) 



{0)t = Tr[Op(t)] 

= ^ e-(^"-^"')*(n5|0|nV)(«V|po|nff), (21) 

nn'gg' 

where \ng) are the eigenstates of H with energies i5„, 
and g labels possible degeneracies. If the long-time limit 
limt_»oo {0)t exists, then it is necessarily equal to the 
long-time average (O), 





= Y,{n9\0\ng'){ng'\po\ng), (22) 

ngg' 



assuming that the limit can be taken termwise (which is 
allowed for the large but finite systems that we have in 
mind). In a steady state the s^^stem is thus described by 
the "diagonal ensemble" Q [II, [H, , 



Pd 



lag 



\ng){ng\po\ng'){7ig'\^Y.P^PoP^- (^3) 



ngg' 



Here P„ — '^g\ng){ng\ is the projector onto the sub- 
space spanned by the eigenvectors corresponding to the 
energy eigenvalue £"„. The statistical operator p^jj^^g cor- 
rectly describes the long-time limit, if it exists, of any 
observable O, i.e., JO) = Tr[Opji^g]. 

The diagonal ensemble correctly yields any stationary 
expectation value, regardless of the transient behavior. 
However, not a lot is gained by regarding Pji^g as a "sta- 
tistical" prediction for the steady state, because each en- 
ergy eigenstate contributes according to the initial con- 
ditions given by po- For a nonintegrable system one 
expects, by contrast, that only a few conserved quan- 
tities such as energy and particle number need to be 
fixed for a successful statistical description in terms of 
thermal Gibbs ensembles, which can emerge from the 
diagonal ensemble by means of "eigenstate thermaliza- 
tion" [li,[3i,[4l|. Similarly one can ask when a GGE ^ 
for an integrable system yields the same prediction as 
the diagonal ensemble . This is discussed in the next 
subsection. 



Gibbs ensemble for an integrable system 



We now consider an integrable system whose time evo- 
lution after the quench is governed by the effective Hamil- 
tonian ([2]). For simplicity we consider two typical cases 
of Hamiltonians only. Either (a) the constants of motion 
Xa have the eigenvalues and 1 and can thus be repre- 
sented by fermions or hard-core bosons, Ta = a^^a^, with 
[Oa, a^]± = Sa/B, (a„)2 = (ajj)^ = 0; or (b) theXa have the 
eigenvalues 0,1,2... and can be represented by bosons, 

= ^i^ai '^itl^ [^a>^^] = ^ap- Examples for case (a) 
are the effective Hamiltonians for hard-core bosons in one 
dimension free fermions with quenched disorder p^ . 
and the 1/r fermionic Hubbard chain (Se cjl|) . whereas 
case (b) apphes to the Luttinger model Q. For both 
cases the Lagrange multipliers Aq, in Eq. ([3]) are then 
given by (a) ln[(X„)o ^ - 1] and (b) ln[(X„)o ^ + 1]. 

The Hamiltonian ^ has the eigenstates \m) with oc- 
cupation numbers Xajm) = ma\m) and energy eigenval- 
ues Em = X]q ^aTna- For simplicity we assume that the 
degeneracy of energy eigenvalues is irrelevant, i.e., the 
observable O or the initial-state density matrix po are 
diagonal in the subspace of eigenvectors \m) with the 
same energy. This assumption will be examined in detail 
at the end of this section. From Eq. the diagonal 
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ensemble and long-time average are then given by 

(O>" = ^(m|0|m)(m|po|m>. (24) 

m 

Is this the steady-state value predicted by the GGE ^1 
We answer this question for two types of observables: for 
case (a) we consider the observable 



A — Aai-'-am t . . . „t 



(25) 



Ql •■ -0:77 

/3i--/3„ 



while for case (b) we allow for powers of the bosonic 
operators and consider (for ri,Sj > 1) 



B 



E 



ibly'---ibLr-ib,j'----{bp^r, (26) 

We assume without loss of generality that -B^^\ . .^'"'J^\. . J" 
vanishes whenever two indices or two indices Pj are 
the same. 

It is straightforward to obtain the long-time aver- 
age (|24l) and GGE average ([3]) of the observables A and 
B by using the occupation number basis |m) and the 
fixed GGE averages (Iq)q = (Tq)o- case (a) we find 



(27a) 
(27b) 



(A) = ™ ( n ) , 



Ql ■ --am 

where we used the identity 



i—l Iq i—1 i—1 



in the second line. In case (b) we have 

OLl ■ ■■CH.-i-n 

= E ^z:: 

Ql ■ ■■Q.rn 
Ql ■ --arn 

= E 



(28) 



UiblY'ibc.Y 

m Vi — l 

n n (la. - k) 

i=l k=0 Iq 



1 = 1 IQ 



(29a) 



n ^.!((i«.)or 



where we used the bosonic operator identity 

(6Lr(^o.r = 'n'(ia. -fc) 

k=0 



(29b) 



(30) 



in the first line, and the identity 



mbiJ'KJ^) = Y{{(blY'ibo.J^)o 

i=l /g i=l 



fc) 



n ( n (la. 

z^l \ fc^O 
m 

n ''.!((ia.)or 



n n!((ia.)Gr 



(31) 



in the second line. Furthermore we defined the permu- 
tation-averaged matrix elements Aa-^...^^ = Tip ("pl)^ 

^QPi - ap„ "qi - q„ ~ Z^P -"api - apm^i-pi - rpm • 

From these results we obtain rather general suffi- 
cient conditions for the validity of the GGE predictions, 
namely the factorization of initial-state expectation val- 
ues of (a) products or (b) polynomials of the constants 
of motion Xq, as follows: 



If ( n ) = n (^'oJo then {A) = . (32a) 

\ 1=1 /o i=l 
I rn ri~\ \ ni - 

If ( n n i^c.-k)) = n ^.!((^a.>or 

\ i=l k=0 /q i=l 

then JB) = {B)^. (32b) 

The ai and in ([32a| and (j32bp are those for which 
^Qi - a,„ and , respectively, are nonzero. The 



criteria ([32|) are the central result of this section, and we 
now discuss their implications in detail. 

First we note that for simple observables, which in- 
volve at most one factor la , the factorizations ([32]) occur 
trivially. This is the case, e.g., for the double occupancy 
in the 1/r Hubbard model, explaining why the GGE (fTOl) 
works in this case. Typical observables of an interacting 
integrable system, however, are often rather complicated 
when expressed in terms of the constants of motion Xa- 
Thus all correlations between the constants of motion 
must vanish in the initial state in order for (j32p to be 
fulfilled. This seemingly restrictive condition can nev- 
ertheless be met, because often the initial state is not 
strongly correlated in terms of the Xa- Moreover, some 
correlations among the Xa are allowed provided that their 
contribution to the sum (P7)) or (^5]) is negligible, e.g., if 
it vanishes in the thermodynamic limit. 

For example, one-dimensional hard-core bosons are 
represented by a free-fermion Hamiltonian. For an alter- 
nating potential (as studied in Refs. @) the initial state 
contains correlations only between the fcrmionic momen- 
tum number operators and Uk+Ti- One can then show 
from ((27l) that the GGE ([3]) makes correct predictions, 
up to finite-size corrections, for observables that are re- 
stricted to a finite region of real space. Another example 
is the fermionic Luttinger model, which maps to a free- 
boson Hamiltonian. For an interaction quench (studied 
in Rcf. |§) the initial state is a product state with corre- 
lations only between the bosonic momentum occupation 
Uq and n_q. 
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On the other hand, correlations between constants of 
motion in the initial state, which remain for all times, 
cannot be described by the GGE ^ , as noted in Ref. [isl . 
However, in interacting integrable systems, such observ- 
ables usually correspond to complicated many-particle 
operators in terms of the original microscopic degrees of 
freedom, and thus arc not measurable in practice. As 
mentioned in the Introduction, the microcanonical Gibbs 
ensemble faces the same problem: one can always con- 
struct fine-grained observables that do not thermalize. 

Can a GGE be improved if it does not yield the correct 
long-time average? The minimal necessary extension of 
the ensemble depends on the observable in question, but 
as one sees from Eq. (f32|) , it always suffices to fix not only 
the expectation values of the constants of motion Xq but 
also the expectation values of all of their products, i.e., 
by using 



Pq oc exp 



(33) 



where the Lagrange multipliers are chosen to fix all prod- 
ucts, Tr[Xaj • • - la^JjQ] = Tr[XQj • • -Xa^Po]- Then it fol- 
lows immediately that Tr[ApQ] JA) and Tr[BpQ] = Jb) 
for the observables considered above. Thus any steady 
state can he described by a sufficiently extended GGE, 
i.e., by fixing sufficiently many products of constants of 
motion. 

While fixing all products as in Eq. ([33|) yields the exact 
long-time average, this extension of the GGE can hardly 
be regarded as a statistical description of the steady state 
(as noted in the previous subsection for the diagonal en- 
semble) , because it uses almost the full information about 
the initial state. In fact, any nondegenerate Hamiltonian 
H acting on a Hilbert space of dimension h has h — 1 
pairwise commuting and linearly independent constants 
of motion fixirig all of them in a GGE recovers the 
diagonal ensemble [40[. For a nonintegrable system one 
can choose, e.g., the projectors onto the eigenstates of 
H [iOl , or linearly independent integer powers of H [ll[ . 
In practice, however, these extensions of the GGE are as 
hard to calculate as the long-time average. 



Degenerate energy levels 



In the previous subsection we assumed that the de- 
generacy of energy levels is irrelevant [as defined above 
Eq. ([23])], which allowed us to move from Eqs. (|22j|23p 
to p4)) . Below we provide an example for which this as- 
sumption docs not hold. In that case the expression (|24p 
for the long-time average cannot be used, and thus nei- 
ther (|27ap nor (|29ap are available. 

We consider a quench to J7 = in a general fermionic 
Hubbard model. Fermions (with spin cr = "f, |) on a Bra- 
vais lattice (with L lattice sites) are prepared in a corre- 
lated unpolarized initial state po with fixed densities n| 



= ni = n/2. The time evolution is governed by the free 
Hamiltonian 



H 



E 



(34) 



where k labels the crystal momentum (we suppress the 
vector notation); periodic boundary conditions are as- 
sumed for simplicity. This Hamiltonian of the form ([2|). 
with the number operators ni;^ = c\^Cf,^ playing the role 
of the constants of motion ■ We are interested in the 
steady-state expectation value of the double occupation 

7i,, 1-714. 

Assuming again that the degeneracy of energy levels 
is irrelevant, we obtain for the long-time average (|24p. 
using the basis \m) = Y{kc,{c\^)'^'"'\0) , 

("*T"u) = E^"^l^ol"^) 72 E ™feT"^fe'i 



kk' 



= Tr 



r 2 



n^nl = — . (35) 



The same value is obtained from the canonical and grand- 
canonical ensemble, and also from the generalized Gibbs 
ensemble which uses the number operators rika as con- 
stants of motion: 



VM)g = ;^E^'^fcT)G("fc'i)( 



kk' 



= I2 E^'^'^-T)o("fe'i)o = "T"i = ^ ■ (36) 

kk' 

Thus we conclude that the double occupation thermalizes 
to the value ri^ /A after a quench to U — in any Hubbard 
model, provided that the degeneracy of energy levels is 
indeed irrelevant. 

Interestingly, this statement disagrees with our exact 
results for the l/r Hubbard chain from Sec. |TT1 when 
quenching from U = cx) to we obtained the long-time 
limit as doo =n'^{3~ 2n)/6 [Eq. US])]. This differs from 
the long-time average ([55)) because the degeneracy of en- 
ergy levels is in fact relevant for this quench: the initial- 
state density matrix does not factorize in the free-fermion 
basis and the linear dispersion Cfc = tk leads to a massive 
degeneracy for the free-fermion energy eigenstates \m). 
Therefore the long-time limit is not given by Eq. and 
does not equal 7i^/4. 

Furthermore this example shows that GGEs based on 
different representations of the constants of motion can 
yield different results. The free-fermion GGE ([55)) pre- 
dicts the wrong value 7i^/4, whereas the GGE ([TO)) which 
uses the effective bosonic representation gives the cor- 
rect value 71^(3 — 2n)/6. In fact the choice of constants 
of motion for the construction of a GGE ^ is always 
ambiguous, as discussed in the Introduction. Neverthe- 
less it is possible to determine the correct GGE a priori. 
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i.e., without knowing the real-time dynamics, by verify- 
ing that the degeneracy of energy levels is irrelevant and 
that the conditions ([32| are fulfilled. 



quantum systems that the GGE prediction equals the 
long-time average, provided that the obscrvablcs or ini- 
tial states arc sufficiently uncorrelated in terms of the 
constants of motion. 



IV. CONCLUSION 

The exact real-time dynamics of the double occupa- 
tion in the fcrmionic 1/r Hubbard chain shows that the 
presence of a Mott gap docs not inhibit the relaxation 
after an interaction quench. Its steady-state properties 
are correctly predicted by generalized Gibbs ensembles. 
Furthermore we showed for a general class of integrable 
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